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ABSTRACT 

Young binaries within dense molecular clouds are subject to dynamical friction from 
ambient gas. Consequently, their orbits decay, with both the separation and period 
decreasing in time. A simple analytic expression is derived for this braking torque. The 
derivation utilizes the fact that each binary acts as a quadrupolar source of acoustic 
waves. The acoustic disturbance has the morphology of a two-armed spiral and carries 
off angular momentum. From the expression for the braking torque, the binary orbital 
evolution is also determined analytically. This type of merger may help explain the 
origin of high-mass stars. If infrared dark clouds, with peak densities up to 10 7 cm -3 , 
contain low-mass binaries, those with separations less than 100 AU merge within about 
10 5 yr. During the last few thousand years of the process, the rate of mechanical energy 
deposition in the gas exceeds the stars' radiative luminosity. Successive mergers may 
lead to the massive star formation believed to occur in these clouds. 

Subject headings: binaries: general — stars: formation — stars: early-type — ISM: 
clouds 

1. Introduction 

The youngest binaries are still embedded in molecular cloud gas. During this phase, the 
orbiting components experience dynamical friction with the surrounding medium. The braking 
torque causes the stars to spiral inward, so that their separation and period diminish with time. 
Such orbital decay may be significant in the very densest molecular clumps. In particular, infrared 
dark clouds, currentl y believed to be the birthsites of massive s tars, have peak number densities 



as high as 10 cm 3 ( Rathborne et al.1 12000 : iBeuther et al.1 120071 ). This impressive fig ure matches 



that in hot molecular cores, already known to contain luminous, high-mass objects (jKurtz et al 



2000). Binary mergers could be occurring in these environments, and could play a role in forming 



the massive stars. 
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There is still no direct evidence that either infrared dark clouds or hot molecular cores contain 
low-mass stars. H owever, observation s show that, with few exceptions, massive stars originate in 
populous clusters (|de Wit et al.ll2005l ). Indeed, one account of h igh-mass star formation is that it 



proce eds through the coalescence of lower-mass cluster members (jBonnell et al.lll99a : IStahler et al 



2000). At even the highest observed cluster densities, collisions between stars, either bare or with 
extended disks, occur at too low a frequency. However, the presence of dense, background gas 
should enhance the capture rate (jBally Zinneckeri 120051 ) . Specifically, the braking torque from 
dynamical friction might not only shrink initially wide binaries, but cause the component stars to 
merge. An accelerating succession of such mergers could then build up massive stars. 

To assess this possibility, the present study calculates the torque on a binary embedded within 
an extensive gas cloud. Although the torque is created by dy namical frict i on, th e usual analysis 
of that effect is difficult to apply in this context. Following iDokuchaevi (119641 ). theorists have 
dete rmined the frictional drag on a single mass traveling th r ough ambient g as in a straight line (see 



also iRuderman Spiegel 



dynamical case (jChandrasekhar 



1971 



Rephaeli fc Salpeteri Il980l : lOstrikerl Il999l ) . Just as in the stellar 
19421). the moving mass draws external matter into a trailing 



wake. It is the gravitational tug from this wake that provides the effective drag. The same essential 
mechanism is at work in the case of a binary. However, the trajectories of both perturbing masses 
are now closed orbits. Further complicating the analysis is the fact that each star draws in matter 
from the wake of its companion. Naive application of the standard dynamical friction formulae 
would grossly misestimate the torque. 

The new approach introduced here concentrates on the fact that the binary as a whole must 
shed angular momentum to the external medium. The actual mechanism is that the orbiting 
stars create an oscillating gravitational potential that torques nearby gas. This fluctuating torque 
generates outgoing acoustic waves, which transport angular momentum. By calculating the total 
angular momentum efflux from the binary, the braking torque may be obtained without considering 
the complex star-gas interaction close to the stars themselves. 

Section 2 below formulates the problem mathematically and derives the governing wave equa- 
tion. The quadrupolar driving potential created by the rotating stars is established in Section 3, 
while Section 4 discusses the physical character of the generated waves. Section 5 presents the 
central result of this investigation. Here the angular momentum transported by the spiral wave is 
obtained, and thereby the torque (see eq. (39)). It is also shown that the wave carries off all the 
mechanical energy released by the shrinking binary. The resulting temporal evolution of the sys- 
tem is considered in Section 6, along with the efficacy of braking in infrared dark clouds. Finally, 
Section 7 compares the results obtained here with the traditional analysis and discusses future 
extensions of this study. 
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2. Formulation of the Problem 
2.1. Length Scales 

We wish to treat the binary as a perturbing mass embedded in otherwise static, uniform gas. If 
the binary itself recently formed, then gas surrounding it would not be perfectly quiescent. Hence, 
our assumed background is a highly idealized representation of a real cloud, or at least that portion 
of a cloud in which we can accurately follow the propagation of acoustic waves from the stars. One 
obvious stipulation is that this region cannot be so close to the binary that cloud gas is infalling 
onto the stars. This requirement sets an inner radius of validity for the analysis, which we take to 
be the sonic point in the Bondi accretion problem: 

r in = GM tot /2c 2 s (la) 
= W0AV(^-) (-^—rY 2 • (lb) 



1 M Q J \2 km s- 1 



Here, M to t is the total binary mass and c s the sound speed of the surrounding gas, assumed to 
be isothermal. In the numerical ev aluation of r\ n , we hay e used for c s the typical observed velocity 
dispersion of infrared dark clouds ( Sridharan et al. 2005l )f*l Inside r- m , where the binary itself may 



reside, gas falls freely onto the central stars. 

If the binary separation, to be denoted as a tot , exceeds r; n , then each orbiting star is sur- 
rounded by its own zone of infall. In that case, atot becomes the appropriate inner boundary. This 
comparison of atot relative to n n effectively distinguishes two cases: "hard" binaries, in which the 
relative speed of the component stars is supersonic with respect to the gas (a to t < Hn) and "soft" 
binaries, for which this speed is subsonic (atot > Hn)- We shall concentrate on the first case, since 
it yields, as we shall see in Section 4, a simplified form of the acoustic disturbance. 

We also require that our propagation region is not so spatially extended that the dominant 
gravitational force is from interior gas rather than the stars. Thus, we are confined within another 
radius r gas , where the mass of background gas rivals that of the stars. For an H2 number density of 
10 7 cm -3 , representating the most compact regions of infrared dark clouds, we find that 1 Mq is 
cont ained in r gas = 1600 AU. Beyond even this point, theory indicates, and observations confirm 



(e.g. iButler Tanl 120091 ) . that the mass density p falls off in response to the self-gravity of the 



cloud as a whole. This outermost radius is the Jeans length Aj, where 

Aj - -£= (2a) 
VG p 



1 We do not know the physical origin of the velocity dispersion, which may represent magnetically mediated 
turbulence, as is believed to be the case in local clouds. We follow the conventional route of assuming that this 
motion is isotropic and contributes an effective pressure acting o n the cloud as a whole. In more complete studies of 
molecular clouds, composite equations of state are emp loved fe.g.lMcKee fc Hollimanlll999h . For the present work, a 
simpler, isothermal equation of state suffices (see also lHuff fc Stahlerll2007 ). 
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9 x 10 3 AU 



n 



10 7 



cm 



-1/2 



2 km s 



(2b) 



Our wave analysis is thus valid within a somewhat restricted, but well-defined, dynamic range. 

By the same token, we cannot model, using this method, binaries embedded within the dense 
cores characterizing low-mass star formation environments (see lHaisch et al.ll2004l . for observations 
of candidate stellar pairs). Here, the entire cloud mass is comparable to the stellar Mtot, and 
fltot < Hri ~ ?~gas ~ Aj. Physically, the dynamics inside such a dense core is dominated by infall 
onto the binary, unless the cloud is stabilized by tension in its internal magnetic field. Acoustic 
waves may propagate outside the cloud, but their source would be both the stars and the dense 
core gas. 



2.2. Perturbation Analysis 

We now introduce our central object, the binary. Let the component stars have masses M\ 
and M2, whose values are in the ratio q, where q < 1. We assume, for simplicity, that the two are 
on circular orbits about their mutual center of mass, with radii a\ and 02- For a given Mtot, q, and 
a tot = a i + a 2i the fixed separation between the stars, we want to ascertain the perturbing effect 
of the binary on the surrounding cloud, taken to have a background density of po- 

In the linear approximation valid for small disturbances, momentum conservation reads 



du\ c 



"-Vpi-V^. (3) 

dt po 

Here we have denoted the small induced velocity by u\, and have used the same subscript for 
perturbations in the density and potential. The superscript on emphasizes that this perturbation 
is due to the stars alone. The equation of mass continuity is, again to linear order, 

^i + p V- Ul = 0. (4) 

In practice, the term V<I>^ in the momentum equation is relatively small. However, this term cannot 
be discarded, as it is the stellar gravity that ultimately drives the outoging waves. 

The next step is to manipulate equations (3) and (4) to obtain a wave equation for p\. Taking 
the divergence of equation (3) yields 

dui = _£i v 2 _ y2$ * 

dt po P 1 
But the time derivative of equation (4) may be cast as 

V — - 1 d 2 Pi 
dt po dt 2 
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Combining the last two equations yields 



$-8*- = (5) 



where p\ obeys Poisson's equation: 



V 2 $^ = 4ttG P 1 . (6) 



Equation (5) is the inhomogeneous wave equation used in all analyses of dynamical friction 
in gases. In previous studies, however, the source density p\ tracks the straight-line trajectory of 
the gravitating mass. After solving the wave equation for p\, the retarding force on that mass is 
calculated. In our case, the gravitating mass is spatially confined. It is thus more convenient to 
idealize p\ as being non-zero only at the origin. This singular, "equivalent density" for the binary is 
derived in the next section, where we also show how its temporal change generates acoustic waves. 



Binaries as Acoustic Sources 



3.1. Potential of the Binary 

We have written the source term of the wave equation (5) in terms of the equivalent density 
p\. To derive this density, we begin with <&l, the potential due to the binary. For any spatially 
co mpact source, the potential is conveniently expressed as a sum over multipoles. Equation (4.2) 
of I Jackson! (|19620 gives the electrostatic expression; this is readily altered to the gravitational one: 

^) = -^G±J2J^^§fi. (7) 



1=0 m=-l 



21 + 1 



Here, the multipole moments are found from the source density (| Jackson! Il962l . eq. (4.3)): 



Qh 



ee Y^O'A') (rj p(x')d 3 x' . 



(8) 



In these two expressions, primed coordinates refer to source points, and unprimed ones to the field 
points where $J is measured. 

Figure 1 shows the binary lying in the x-y plane of our coordinate system, with its center of 
mass situated at the origin. Also indicated are two of the three coordinates (r, 8, 4>) of a field point. 
At time t, mass M\ is a distance r' = a\ from the origin, and has rotated from the x-axis by the 
azimuthal angle 4>' = tot, where oj is the binary's angular speed. The other mass M2 is located at 
r' = 02 and eft' = tut + it. Both masses have polar angles 9' = tt/2. 



Evaluating the qi m from equation (8), we find that the monopole (/ = 0) term is 

<?oo = Yn, 
M tot 



■00 Mt + r^(|,wt + 7r) M 2 



4tt 



0) 
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In contrast, all dipole (I = 1) terms vanish. For example, 



9ii 




— (a\M\ — CL2M2) exp (— iut) 



= 0, (10) 

by definition of the center of mass. 

The quadrupole (I = 2) terms are the most interesting for our purpose. For m = 0, we have 



920 = 4vi J ' 



(11) 



where I = Mi a? + M2 a| is the binary's moment of inertia about its axis of rotation. For m = 1, 



while for m = 2, 



921 = , 



1 /15 T 
922 = ^ Y — I exp (-2iut) 



(12) 
(13) 



To complete the set, equation (4.9) of I Jackson! ([1962) tells us that 



92,-1 



"921 



(14) 



and 



92,-2 — 922 

1 /15 
4V 2^ 



I exp (+2«c>i) . 



(15) 



Substitution of these coefficients into equation (7) yields the potential through quadrupole 
order: 



which simplifies to 



GM tot GI /3 2 . 1 

5- - cos 6 — - 

r 2r 3 \2 2 

3GI 

~8^~ 



sin 2 6> {exp [+2i(o;t - <j>)\ + exp [-2i(wi - </>)]} 



GMtot _ G/ 
t 2 r 3 



3 

- sin 2 9 cos 2(wi 



-cos 



(16) 



Octupole and higher-order terms also exist, but they are of smaller magnitude. 
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3.2. Equivalent Density 

We now investigate what distribution of matter generates the potential <I>|. This equivalent 
density p* is a simplified representation of the binary that gives rise, through its fluctuating po- 
tential, to the same acoustic waves. Direct calculation reveals that V 2 = outside the origin. 
Thus, the equivalent density is confined to r = 0. We already know that the monopole term in <&l 
is generated by a mass M tot at the origin. We shall denote this portion of the equivalent density 
as pl , and caution the reader not to confuse the double subscript with that in the quadrupole 
moments q\ m . We then have 

Plo = M tot S(x')S(y')5(z') . (17) 

Higher multipole terms in the potential similarly have sources, i.e., further contributions to the 
equivalent density that are combinations of 5 functions and their derivatives, as discussed in 



Morse k Feshbachl (|1953l . p. 1278). 



The quadrupole portion of the equivalent density is illustrated in Figure 2. Here we show a 
symmetric placement of point masses that generates the corresponding part of $*. We replace the 
real binary by two equal masses lying along the binary's azimuthal direction within the x-y plane. 
Each point is at the same small distance e from the z-axis, and the value of each individual mass 
is I/2e 2 . In addition, we place a third point source, with negative mass —I/e 2 , at the origin itself. 

This configuration clearly has no monopole moment (?oo, since its total mass vanishes. It may 
also be verified that all dipole moments are zero. The quadrupole moment for m = is 

in agreement with equation (11). The other quadrupole moments similarly match. 

Referring to Figure 2, the corresponding portion of the equivalent density, which we denote as 



Pl2. 1S 



P12 = 7T" o [Hx' — e cos Lot) 5 ( y' — e sin ujt) 5{z' 



I 

2? 

+ 5(x + ecoswt) 5(y' + esinwi) 8(z') 

-2 5(x')5(y')5(z')] . (19) 



We now expand the arguments of the 5 functions to second order in e. For example, we write 

S(x' — ecosojt) = S(x') — ecosojt d'(x') + — cos 2 ujt 5"(x') . 

We substitute these series back into equation (19) and multiply them out, retaining only terms 
through order e 2 . Many terms cancel, leaving 

Pl2 = P*A + Pb + Pc . ( 2 °) 
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where 



p* A = J - cos 2 ut 5"{x) 5{y) 5{z) (21a) 
Pb = J - sin 2 u;t S(x) 8"(y) 5(z) (21b) 
p* c = - sin2wi 5'(x) S'(y) 5(z) . (21c) 



Notice that p\ 2 varies temporally with a frequency twice that of the binary itself. The reason 
for this doubling is that the simplified configuration in Figure 2 repeats itself every half period. 
That is, any difference in the component masses manifests itself only in higher-order multipole 
terms. 

The full equivalent density is 

Pi = Pxq + P*12 • (22) 
Since we have purposefully designed p\ to have the same monopole and quadrupole moments as the 
real binary, the generated potential must be given by equation (16). Nevertheless, an important 
check is to verify this fact directly, through integraton of Poisson's equation. We carry out the 
relevant calculation in the Appendix. 



4. Character of the Wave 

4.1. Prologue: Static Solution 

To find pi, we substitute our expression for p\ into the wave equation (5). The static portion 
of p*, given by p\$ in equation (17), similarly generates a static contribution to pi, which we denote 
as p\. That is, p\ obeys 

2 AirGpoMtot . . . . , . 
V p x = 5(x) 5{y) 5{z) . (23) 

c% 

Because of the linearity of the full wave equation, p\ can (and should) be added to the time-varying 
part we shall derive shortly. 

Equation (23) is Poisson's equation, describing the potential from a point mass located at the 
origin. The value of this fictitious "mass" is — poM to t/c 2 . We may immediately write down the 
solution: 

= po G M tot 
Cg r 

A similar, but modified, exp ression applies t o the wake created by a moving mass in the traditional 



dynamical friction problem. lOstrikerl ()1999l ) has noted that equation (24) is the linear approxima- 



tion to the full density in a hydrostatic envelope surrounding a gravitating mass M tot : 

p{r) = po exp \-^r ) • ( 25 ) 
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The sharp density rise interior to r = r- m predicted by equation (25) does not actually occur. As 
we have described, gas in this region instead goes into free-fall collapse onto the stars. 

The linear result given by equation (24) is thus the appropriate form for the static density 
enhancement El However, p\ makes no contribution to the angular momentum transport. In calcu- 
lating the latter, we shall be multiplying the density by the induced, azimuthal velocity. Since the 
latter oscillates sinusoidally, the product vanishes over a period. We therefore turn to the oscillating 
density perturbation. 



4.2. Wave Density 

From now on, we may omit p\ Q when considering the equivalent density. The time-varying 
density perturbation, which we shall continue to denote simply as p\, obeys 

Vpl ~cl^ = c^ Pl2 ' (26) 

where p\ 2 is given by equations (20)-(21). We proceed by finding those parts of p\ (denoted pa, 
etc.) generated by each additive component of p\ 2 - Linearity of the wave equation ensures that we 
can add these individual solutions to obtain the full one. 

Consider first the functions T>a, T^b, and T>c obeying 

v v A — ^ ~p^r = 2 — cos ut °y x > °(y> °( z ) ( 27a ) 

c s at c s 
V ^"^^ = - 27rG / 0/ S m 2 ^ 5(x)S(y) 5{z) (27b) 

C% Ot Z Cg 

V2V c--2^ = - 2 -^^^t5{x)5{y)5{z). (27c) 

Cg Ot z Cg 

If we can find these three functions, then differentiation of their governing wave equations reveals 
that 

d2VA /no \ 

PA = -q^T ( 28a ) 



PB = ^ (28b) 



PG = (28C) 



2 Actually, equation (24) only gives the static density perturbation to leading order. As we shall see in Section 4.3, 
another term that dies off as r -3 is needed to balance the static part of the quadrupole potential. 
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Each of the wave equations (27a)-( 27c) may be s olved using the retarded Greens function. 
Applying equations (6.54) and (6.66) of I Jackson! (|1962l ). and integrating over the 5 functions, we 
find 



(29a) 
(29b) 
(29c) 



V A = 


G Po I 

2c 2 s r 


cos 2 (ut — 


kr) 


v B = 


Gp Q I 
2 c 2 r 


sin 2 (ut — 


kr) 


V c = 


G Po I 

2 c 2 r 


sin 2 (ut — 


kr) 



where the wave number k = u/c s . 

In taking spatial derivatives of these last three expressions, we utilize the fact that we are in 
the far-field limit (kr S> 1). To see this, note first that 



k 2 n 2 
K "tot 



G Mi 



tot 



c 2 a to t 



2 r in 

"tot 

Since our field point is located well outside r- m , we have, for hard binaries, 

k r k r-. u ^> k fl+ n + . 



We conclude that 



; 2 2 ^ ^ Tin 

k r » 

"tot 



> 1 . 



Following the usual practice in acoustics (e.g. lLighthilllll978l . Chapter 1), we apply spatial derivatives 
only to the phase (ujt — kr). Derivatives of the prefactors of T>a, etc. involving r are smaller by 
one or two powers of (kr)^ 1 . 

Differentiation of equations (29a)- (29c), under the far- field approximation, yields expressions 
for the density components: 



PA 



PB 



PC 



— ^ sin 2 # cos 2 (/> cos 2 (cut — kr) 



+ 



Gpplk 2 

c 2 r 
Gpplk 2 



sin 2 # sin 2 cos 2 (ut — kr) 



c? r 



sin 9 sin 2(f) sin 2 (ut — kr) . 



Adding these gives the full density perturbation: 



Gpplk 2 . 2 

pi = k sin cos 2(ut — kr 

cjr 



(30a) 
(30b) 
(30c) 

(31) 
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lt is important to understand, in a qualitative sense, the amplitude of the density perturbation 
in equation (31). The relative perturbation, pi/po, created by a simple point mass (a monopole) 
is of order r; n /r, according to equation (24). However, our oscillating density perturbation is 
quadrupolar. Thus, the monopole result must be multiplied by two powers of k a to t- The amplitude 
in equation (31) is indeed of order (r in /r) (fca to t) 2 - 

As expected, the perturbation is an acoustic wave that travels radially outward with phase 
velocity u/k = c s . At any time, the phase of the wave is also dependent on <j>. In fact, equation (31) 
reveals that the disturbance may also be viewed as a trailing, two-armed spiral wave, with a 
latitude-dependent amplitude that peaks at the equator (6 = tt/2). Since fcr> 1, the spiral is 
tightly wrapped, with a relatively small pitch angle. 

Figure 3 illustrates the basic geometry of the wave. Shown are wavefronts (surfaces of constant 
phase) for the two spiral arms in the equatorial plane. If we trace one arm around the circle, the 
radius of the front increases by A = 2ir/k. However, because a second arm is interleaved, the 
actual radial wavelength of the disturbance is A/2, with an associated wavenumber of 2 k. The 
perturbation's angular frequency is 2u, so the outward velocity is again 2uj/2k = c s . 



4.3. Induced Velocity 

We next determine the velocity created in the gas by the passing wave. Taking the curl of the 
momentum equation (3), we find that 

|(Vx« 1 ) = 0. 

Thus, the induced vorticity is independent of time, and is zero for oscillatory motion. It follows 
that the velocity may be written as 

«! = Wi , (32) 
where ipi is the velocity potential. From the mass continuity equation (4), ipi obeys 

V^ = -i%. (33) 
Po at 

If we assume that ipi depends on the same phase as p\, then the dominant contribution to 
V 2 ^i in the far-field limit is simply —Ak 2 ip\. Using pi(t) from equation (31), we find that 

V>i = — ^ sin 2 sin2(u;t - kr - <f>) . (34) 
2 cj r 

Finally, we may read off from equation (32) the velocity components: 

G I J 1 

u r = ^ — sin 2 # cos2(wt — kr — <f>) (35a) 
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GIlo 
2 c 2 . r 2 
GIlo 



sinO cos2(o;t — kr — <j)) . 



sin 2$ sin2(wt — kr — <j)) 



(35b) 



(35c) 




In deriving u r , we again applied the radial derivative to the phase only. We see also that 
this velocity component is larger than the other two by a factor of order ur/c s = kr 3> 1. Such 
dominance of the radial velocity is expected for a wavefront with small pitch angle. While relatively 
small, the ^-component is critical for angular momentum transport. 

We may also obtain u\ directly from the momentum equation (3). The gradient of the 
monopole contribution to Q\ is balanced by the static density perturbation previously derived. 
Equation (16) shows that there is also a static part of the quadrupolar potential. This is balanced 
by a smaller term in the static density perturbation. That is, the total static perturbation, to 
quadrupole order, is 



The remaining, oscillatory part of $^ contributes in principle to the fluctuating velocity. How- 
ever, if we actually compare its gradient to the force associated with the pressure gradient, we 
find the (^-gradient to be smaller by several powers of (kr)^ 1 . In the far-field limit, therefore, 
the velocity is actually generated only by p±, as equation (33) already indicates. If we express the 
oscillating part of p\ as the real part of a complex exponential, and set 



then we may solve the momentum equation for u\ itself, obtaining the same velocity components 
as above. 



Imagine surrounding the binary with a spherical shell of radius r. We wish to determine the 
outflow of angular momentum through this shell. For the acoustic wave to have the properties we 
ascribed to it, r must lie between r- m and r gas . The z-component of specific angular momentum 
at any point on the shell is rsvuOus, with ua, given by equation (35c). Additionally, the mass 
flux through this same point is pu r , where u r is taken from equation (35a). The flux of angular 
momentum, which we denote as j, is therefore 




(36) 



5. Angular Momentum and Energy Transport 



5.1. Braking Torque 



3 = 




(37a) 



sm 



H cos 2 2(ut - kr - (f>) . 



(37b) 
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In the last expression, we have used for p its equilibrium value pq. Since the product rsm8u r u^ 
already falls off as r~ 2 , any density variation that declines with radius does not appear in the total 
angular momentum efflux, integrated over the sphere. 

By angular momentum conservation, this outflow, which we denote J, must also be —T, where 
r is the torque exerted on the binary by surrounding gas. That is 

r = -r 2 / d(9sin6> / d(f> j . (38) 



o jo 



Using 

r>27T 



/ dO s'm 5 e / d(j) cos 2 2(ut - kr 
Jo Jo 



16 7T 



we arrive at our main result: 

„ 16 7T Ll> 3 



A striking aspect of the torque is its high sensitivity to the sound speed c s . It is more difficult 
to gather hotter gas into the wakes that actually provide the gravitational tug on the orbiting stars. 
An inverse dependence on c s is also present in the expressions for p\ and U{ (see equations (31) and 
(35a)-(35c)). In any event, the sensitivity of T to c s means that quantitative conclusions regarding 
astrophysical effects of the torque are necessarily rather imprecise. 



5.2. Energy Loss 



The outgoing acoustic wave transports not only angular momentum, but also mechanical 
energy. We first note, from equations (35a)- (35c) and the succeeding comments, that the kinetic 
energy density in the far field simplifies to 

\^p(u 2 r + uj + uj) -» ^p u 2 . 

Now the total energy density of any acoustic wave, including t he component a ssociated with com- 
pression by the enhanced pressure, is twice this kinetic value (|Lighthilllll978l . Section 1.3). Since 
the wave travels radially outward at the sound speed, the total energy flux past any point is 



£ 



Pq U r C s 



(40) 



It is instructive to compare this result with j, the angular momentum flux in equation (37a). 
From equations (35a) and (35c), we have 

U <P = I • 

ijj rsm0 

Thus, the angular momentum flux may be written as 



J = —pou r , 



(41) 
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where we have again replaced the density by its equilibrium value. Integrating £ and j over the 
entire shell, we obtain a relationship between the global energy loss rate E and J: 

E = ujj . (42) 

Since both the energy and angular momentum are being extracted from the binary, the same 
relationship between their loss rates should apply to that system. We now show that this is the 
case. Referring back to Figure 1, the angular momentum of the binary is 

Jbin = Mi a\ uj + M2 0,2 oj 

= Iuj . (43) 

1 T 2 G Mi M 2 



The binary's total energy is 



But we also have 



2 atot 

= ^ . (44) 



"tot 

From this equation, applying standard manipulations, we find that the potential energy is 

G Mi M 2 o 
= — 1 uj , 

«tot 

so that 

E hhl = ~Iuj 2 . (45) 
Comparison with equation (43) reveals that 

^bin = -|«/bi„- (46) 

We next relate the temporal change of Ju n to that of uj. We first note that the binary's moment 
of inertia may be written in terms of the separation a tot and mass ratio q: 

1 = jrfw Mtot at2ot ' (47) 

Thus, if we again use equation (44) to eliminate a to t, the angular momentum may be written as 

Jun = M t 5 f W "V3 . (48) 

During contraction of the binary, therefore, 

Jbin 3 UJ 
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Taking the time derivative of equation (46) and applying equation (49) now gives 



W T 
2 bin 


UJ 

2 


Jbin 




LO 


Jbin 


T bin ~ 


2 


U Jbin ■ 







(50) 

As claimed earlier, the energy and angular momentum of the binary change at the same relative 
rates as these same quantities in the outgoing wave. 

The rate of energy transport by the wave can be recast in another way that provides a check 
on our derivation. We first note, after applying equation (47) to the negative of equation (39), that 



J 



16 vr 



Po (GM tot ) 2 a 



tot 



15 (1 + q) 4 c 5 s 

If we then use equation (44) to eliminate G M tot in favor of atot and lo, we obtain 



(51) 



J 



16 7T 



15 (1 + q) A c 
Thus the energy emission rate can be written as 



Po a 



10 

tot 



E 



16 vr 



q 



15 (1 + q) A c; 



-r Po a 



10 



(52) 



For fixed a to t 5 the binary components' relative speed scales with u>. This last expression thus 
reproduces the fact that the acou stic energy rad iated by a quadrupolar source increases as the 
eighth power of the Mach number (|Lighthilllll952l ). 



6. Binary Evolution 

We are now in a position to follow the binary's orbital decay in time. For this purpose, we use 
for the torque the negative of J in equation (51). However, it is now appropriate to eliminate a tot , 
again employing equation (44): 

We set this torque equal to the temporal derivative of Jbin, as given in equation (48). Rearrangement 
gives an equation for the evolution of uj: 



16 vr q ( G M tot \ 5/3 r, a ,- , . 
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Equation (54) is readily integrated. If ujq is the initial, nondimensional angular rotation rate, 

then 

t \ ' 3/2 



u = Uq[1-—] , (55) 



where the coalesence time t r . is 



15 (1 + g) 2 1 /GMtot\- 5/3 _ 2/3 , 

*c = ™ ~ 3— ^0 56a ) 

32 vr g p G V 4 J 

= 2 X 10 yr (10^=3) (2kn7^J {iM^j [WYr) ■ (56b) 

At time t c , uj diverges and the binary has contracted to zero separation. In our numerical evaluation 
of this time, we have set q = 1 and used the initial binary period Po i n place of the angular velocity 
luo. For M t ot = 1 Mq, a period of 10 3 yr corresponds to a to t = 100 AU@ 

Finally, we may determine the mechanical energy release of the decaying binary as a function 
of time. Setting E = —ujT and taking T from equation (53), we have 



16 » f—^^.au 

C 

= ^o(l-|;)" 2 . (57b) 

Here, we have supplied the time dependence of uj from equation (55). The constant Eq is 

■ 16 vr q 2 4/3 f GM tot \ W/3 5 

* = it jiTtf u ° v^rJ P0Cs ' (58) 

and has the value 9 x 10 -4 Lq for our fiducial parameters. Evidently, the energy release does not 
rival the radiative loss from the stars themselves until very late during the inspiral, when t is within 
a few percent of t c . At this point, the stars are separated by several AU. 



Discussion 



We may compare, at least in a qualitative manner, our derived torque w ith that indicat ed by 
the traditional theory of dynamical friction. According to equation (12) of lOstrikerl (| 19991 ). the 
retarding force on a mass M moving at speed V through a cloud of density po is 



Fdf 



4vr(GM) 2 p 
V~ 2 



1 



3 The evolutionary equation (54) neglects accretion from the external medium. A numerical estimate shows that 
the mass gain is not major for our adopted parameters, but it should be included in a more complete analysis. 
According to equation (56), any increase of Mtot shortens the coalescence time t c . 
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The factor X, essentially a Coulomb logarithm, is a nondimensional function of the Mach number 
V/c s and the time since the mass first entered the cloud in question. Since I will generally be 
of order unity, we may ignore it, along with other such factors, in the dimensional argument that 
follows. 

To apply this formula to the binary problem, we interpret V as the components' relative 
velocity V ve \, which is u atot- Then the torque is of order atot ^df, so that 

1 np ~ = . 

Here, the moment of inertia has been approximated as / ~ M a 2 ot . Comparison to equation (39) 
shows that Tdf is the true T multiplied by a factor (c s /VJ- e i) 5 - This factor can be much smaller 
than unity for the hard binaries of interest. On the other hand, the nondimensional I formally 
diverges for a velocity of c s . The traditional theory is thus unreliable in this context. 

However, it is not difficult to envision circumstances in which the present theory requires 
modification. From equation (35a), the Mach number associated with the radial velocity amplitude 
is 

GIlo 2 



II r 



c, cf r 



^rel A '"in 



C, / r 



(59) 
(60) 



Since r m <C r, the induced velocity is normally subsonic. For binaries that are initially very hard, 
or late during the inspiral of any system, u r throughout the far field becomes supersonic. The 
disturbance then changes character from an acoustic wave to a tightly wound spiral shock. In 
the same regime, the density perturbation p\ is comparable to or even exceeds po, so that a fully 
nonlinear treatment is necessary. 

Additional modification of the theory would be required by the inclusion of a finite eccentricity 
in the binary orbits. While the system would still be periodic, the equivalent quadrupolar source 
would exhibit variation over a continuous distribution of frequencies. The same frequency distribu- 
tion would then appear in the transmitted waves. It would be interesting to recalculate both the 
torque and energy loss under this more general condition and thereby follow the evolution of the 
eccentricity during orbital decay. 

Returning to the theory's main astrophysical application, our result for the coalescence time is 
a reasonable one that adds credence to the underlying picture. The nearest and best-studied region 
of m assive star format ion is the Orion Nebula Cluster, whose general population formed 1-2 Myr 



any precision. 



ago (Hillenbrand 119971). We do not k now the age of the high-mass members, the Trapezium, with 



Palla &: Stahlerl (|200ll ) have argued that they are relatively young, of order 10 5 yr, 
based on the location of BM Ori and its binary companion in the HR diagram. Assuming that 
the Trapezium stars, along with their companions, are coeval, their inferred age provides an upper 
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bound to the formation time scale. We are thus encouraged by this matching of times. We stress, 
however, that our expression for t c in equation (56) varies inversely with the imprecisely known 
ambient density n. 

If massive stars indeed coalesce over such a period, perhaps they do so through an accelerating 
sequence of binary mergers. The binaries themselves might be created and disrupted rapidly out 
of the dense gas, l eading to a statistically stable period distribution, as originally envisioned by 



Lvnden-Belll ()1969l ). Within this picture, one could in principle determine the mass distribution of 



the growing population of coalesced objects, thereby advancing the theory another significant step. 



This project was originally inspired by extensive discussions with Avery Broderick concerning 
the formation of massive stars. Kevin Bundy provided useful comments on a preliminary draft of 
the manuscript. The author was partially supported by NSF grant AST-0908573. 
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A. Gravitational Potential from the Equivalent Density 

We wish to verify that p\ in equation (22) indeed generates the gravitational potential <3?| in 
equation (16). It is evident that the monopolar portion p\ Q in equation (17) does yield the term 
G Mt t/r. We therefore focus on p\ 2 , given in equations (20)-(21), and verify that 

**i2(«) = ~G [ ^\ d 3 x> , (Al) 
J \x — x'\ 



where 



3 /3 1 

-sin 2 # cos2(wi — (/))— ( -cos 2 # — - 



(A2) 



The derivation proceeds in a manner analogous to the calculation of p\ in Section 4.2. We first 
find three functions Ta, Tb, and Tc obeying 

V 2 T A = 2ttGI cos 2 cut 5 {x) 5 (y) 5 (z) (A3a) 
V 2 JT B = 2vrG/ sin 2 ^ S(x)6(y)6(z) (A3b) 
V 2 ^" c = 2n G I sin2 tat 5(x)5(y)5(z) . (A3c) 



If we further define 



^ - ^ (A4a) 



dx 2 
dy 1 



**b = (A4b) 



then the combination 



C dxdy ' ( c ) 



$* 2 = <D^ + $*, + $5 , (A5) 



obeys Poisson's equation with p* 2 as the source density. Equivalently, <&\ 2 is the solution of equa- 
tion (Al). 

The functions Ta, Tb and Tc are all solutions of Poisson's equation with a central point mass: 

G I 

Ta = —— — cos 2 ujt (A6a) 
G I 

T B = sin 2 ut (A6b) 

2r 

G I 

T c = -— sin2o;i . (A6c) 

By successive differentiation of 1/r, we find 

GJ / 1 3x 2 



*A = ( ~f z + ^ ) cos 2 ^ (A7a) 
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$*B = ~ (-\ + %)s\n 2 ut (A7b) 



Adding these components yields 



which simplifies to equation (A2). 



2 V r^ r^ 

*h = -^^sin2u,i. (A7c) 



G I 

$* 2 = -_ [cos 2 ojt (-1 + 3 sin 2 9 cos 2 (f) 

_■ 7^ 

+ sin 2 u^ (-1 + 3sin 2 6» sin 2 (f) 

+ sin2u;t (3sin 2 # sm.<p cos^)] , (A8) 
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Fig. 1. — Basic binary parameters. Two stars, M\ and M2, orbit their common center of mass, 
which is at the origin of the coordinate system. The circular orbits, of radii a\ and 02, respectively, 
lie in the x-y plane. At the time shown, M\ has swept out the angle <f>' = cot from the x-axis. 
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Fig. 2. — Equivalent system of masses that gives the same gravitational potential, through 
quadrupolar order, as the actual binary. Two equal masses lie at a distance e on either side of 
the rotation axis, while a negative mass of twice the magnitude lies at the center. 
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Fig. 3. — Geometry of the acoustic wave. The radial wavelength of a single spiral arm is A, as 
shown. However, because a second arm is interleaved, the actual distance between successive wave 
crests is A/2. The central binary is rotating counterclockwise at angular speed uj. 



